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ABSTRACT 

For dyanmic background models of Goldreich & Weber and Lou & Cao, we examine 
three-dimensional perturbation properties of oscillations and instabilities in a gen- 
eral polytropic homologously collapsing stellar core of a relativistically hot medium 
with a polytropic index 7 = 4/3. Perturbation behaviours, especially internal gravity 
g— modes, depend on the variation of specific entropy in the collapsing core. Among 
possible perturbations, we identify acoustic p— modes and surface f— modes as well as 
internal gravity g"*"— modes and g~— modes. As in stellar oscillations of a static star, 
we define g"*"— and g~— modes by the sign of the Brunt- Vaisala buoyancy frequency 
squared N"^ for a collapsing stellar core. A new criterion for the onset of convective 
instabilities is established for a homologous stellar core collapse. We demonstrate that 
the global energy criterion of Chandrasekhar is insufficient to warrant the stability of 
general polytropic equilibria. We confirm the acoustic p— mode stability of Goldreich 
(fe Weber, even though their p— mode eigenvalues appear in systematic errors. Unstable 
modes include g~— modes and sufficiently high-order g'^— modes, both corresponding 
to convective core instabilities. Such instabilities occur before the stellar core bounce, 
in contrast to instabilities in other models of supernova (SN) explosions. The break- 
down of spherical symmetry happens earlier than expected in numerical simulations 
so far. The formation and motion of the central compact object are speculated to be 
much affected by such g— mode instabilities. By estimates of typical parameters, unsta- 
ble low-order I = 1 g— modes may produce initial kicks of the central compact object. 
Other high-order and high-degree unstable g— modes may shred the nascent neutron 
core into pieces without an eventual compact remnant (e.g. SN1987A). Formation of 
binary pulsars and planets around neutron stars might originate from unstable I ~ 2 
g— modes and high-order high-degree g— modes, respectively. 

Key words: hydrodynamics — instabilities — stars: neutron — stars: oscillations 
(including pulsations) — supernovae: general — waves 



and numerical simulations. As nuclear fuels eventually be- 
come insufficient in the stellar core, the process of core- 
collapse SNe signaling the demise of massive progenitors 
(e.g. red and blue giants) may be conceptually divided into 
three stages of core collapse, rebound shock and neutrino 
heating (e.g. Burrows et al. 1995; Janka & Miiller 1996). 

The fortuitous detection of neutrinos from SN1987A 
(Hirata et al. 1987; Bionta et al. 1987; Koshiba 2009 pri- 
vate communications), bolsters such a scenario framework 
in part or as a whole. Optical observations before SN1987A 
revealed its progenitor as a blue giant star in a mass range 
of ~ 16 - 22M0 (e.g. Arnett et al. 1989). At the time of 
SN1987A explosion, twenty neutrinos in the energy range of 
~ 6 — 39 MeV were intercepted within ~ 12 s, confirming 



1 INTRODUCTION 

Supernovae (SNe), hypernovae and a few detected SNe asso- 
ciated with long gamma-ray bursts (GRBs) serve as impor- 
tant cornerstones of several major branches in astrophysics 
and cosmology. Physical mechanisms and outcomes for such 
violent explosions of massive stars have been actively pur- 
sued for decades. Hydrodynamics and magnetohydrodynam- 
ics (MHD) together with simplifying approximations and in- 
creasingly sophisticated microphysics have been invoked to 
model various aspects of SNe in both analytic treatments 
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the occurrence of neutronization. The timescale of neutrino 
emissions was consistent with the prediction for neutrino 
trapping inside an extremely dense collapsed core. The total 
neutrino flux was consistent with energetic neutrinos carry- 
ing off the binding energy during the core neutronization 
(Chevalier 2009 private communications), even though no 
signals of a neutron star (e.g. a pulsar) or a black hole are 
detected (McCray 2009 private communications). 

In spite of extensive research on analytic and numerical 
studies of SNe over several decades (e.g. Goldreich & We- 
ber 1980 - GW hearafter; Yahil 1983; Bruenn 1985; Bruenn 
1989a, b; Herant et al. 1995; Janka & Miiller 1995, 1996; 
Fryer & Warren 2002, 2004; Blondin, Mezzacappa & De- 
Marino 2003; Blondin & Mezzacappa 2006; Burrows et al. 
2006, 2007a, b; Lou & Wang 2006, 2007; Wang & Lou 2007, 
2008; Lou & Cao 2008; Hu & Lou 2009), a few major issues 
in the SN model development remain to be explored (see 
Burrows et al. 2007 for a recent review) . Among these theo- 
retical challenges, the dynamics of core-collapse stage inside 
the progenitor and the possibility of convective instabilities 
during this phase are the main thrust of this paper. 

The simplest hydrodynamic model to describe a core 
collapse is a one-dimensional radial contraction with spher- 
ical symmetry under the self-gravity. In analytical model 
analyses, approximations to the equation of state (EoS) for 
gas medium are necessarily introduced. It can be shown in 
statistical mechanics (e.g. Huang 1987) that a relativistic hot 
Fermi gas with a temperature much lower than the Fermi 
energjijcan be modelled by a simple 7 = 4/3 polytropic EoS 
with 7 being the polytropic index (i.e. the rest mass of a sin- 
gle particle <C the kinetic energy of a particle <^ the Fermi 
energy) . This approximate EoS also gains support in numer- 
ical simulations (e.g. Bethe et al. 1979; Hillebrandt, Nomoto 
& Wolff 1984; Shen et al. 1998). For instance, Bethe et al. 
(1979) concluded that as the neutrino trapping occurs, rela- 
tivistic electrons, high-energy photons and neutrinos mainly 
contribute to the total pressure within a collapsing stellar 
core under gravity. As nuclei start to "feel" each other at a 
high density reaching up to ~ 2.7 x 10^'*g cm~^, the stiffness 
of nuclear matter adjusts the polytropic index 7 to ~ 2.5. 

There are two kinds of polytropic approximations. One 
is the conventional polytropic EoS with P = np^ {P and p 
are pressure and mass density, respectively) where k remains 
constant in space and time. The other is a general polytropic 
EoS where k remains constant along streamlines, i.e. 



where u is the bulk flow velocity. The former is actually a 
special case of the latter, although we treat them separately. 
It is easy to prove that for a relativistic hot Fermi gas with 
a temperature much lower than the Fermi energy, the spe- 
cific entropy of the material is a function of log(P/p^) (see 
footnote 1). In other words, the former EoS is interpreted 
as a constant specific entropy while the latter bears a vari- 
able distribution of specific entropy. The latter general poly- 

^ The Fermi energy is given by Ep = [3h^ pY^/ (iivmp)]^^^ c = 
30(pi4 Ye)^''^ MeV where is the number of electrons per 
baryon, h is the Planck constant, p is the mass density, mp is 
the proton mass, c is the speed of light, pn is the mass density 
in unit of 10^^ g cm~'^. 



tropic EoS obeys the conservation of specific entropy along 
streamlines. The gas dynamics involves the competition be- 
tween self- gravity and pressure gradient force. Such models 
give rise to various self-similar solutions by which the flow 
system partially loses its memory of initial and boundary 
conditions (e.g. Larson 1969; Pension 1969; Shu 1977; Cheng 
1978; GW; Yahfl 1983; Suto & Silk 1988; Lou & Shen 2004; 
Yu et al. 2006; Lou & Wang 2006, 2007; Wang & Lou 2007, 
2008; Lou & Cao 2008; Hu & Lou 2009). For core-collapse 
SNe, GW studied the homologous stellar core collapse of a 
7 = 4/3 conventional polytropic gas. They concluded that 
starting from a static core, if the pressure is reduced by no 
more than ~ 2.9%, the central core will evolve into a ho- 
mologously collapsing phase. This fraction is much less than 
~ 26% as indicated by Bethe et al. (1979). Yahil (1983) 
extended GW analysis to 7 < 4/3 conventional polytropic 
cases and noted the existence of an outer envelope moving 
inwards with a supersonic speed. One important difference 
between GW and Yahil (1983) is that solutions of the former 
have outer boundaries with zero mass density there while 
those of the latter extend to inflnity (Lou & Cao 2008). 

With a similarity transformation (Fatuzzo et al. 2004), 
Lou & Cao (2008) extended homologous core collapse to gen- 
eral polytropic cases. This is a substantial theoretical devel- 
opment of the model framework because several studies (e.g. 
Bethe et al. 1979; Bruenn 1985, 1989b; Burrows et al. 2006; 
Hillebrandt et al. 1984; Janka & Miiller 1995, 1996; Shen et 
al. 1998) on EoS with microphysics during SNe do suggest 
variable specific entropy depending on physical conditions, 
including density, temperature and metallicity. 

Meanwhile, extensive numerical simulations show that 
spherically symmetric models cannot initiate SN explosions 
with an energy of ~ lO'^^ ergs (e.g. Janka & Miiller 1995, 
1996; Kitaura, Janka & Hillebrandt 2006) partly because 
the SN explosion energy appears insufficient and partly be- 
cause instabilities occur in multi-dimensional simulations. 
The role of instabilities and symmetry breaking in SNe has 
now been emphasized (e.g. Burrows 2000, 2006). Along this 
line, various instabilities were proposed (e.g. Goldreich, Lai 
& Sahrling 1996; Lai 2000; Lai & Goldreich 2000; Murphy, 
Burrows & Heger 2004; Blondin et al. 2003; Blondin & Mez- 
zacappa 2006) and several mechanisms may provide seed 
fiuctuations before and during SN explosions (e.g. Bazan 
& Arnett 1998; Meakin & Arnett 2006, 2007a, b). Prior to 
the onset of a core collapse, the so-called "e-mechanism" 
(e.g. Goldreich et al. 1996; Murphy et al. 2004) may lead to 
g— mode overstabilities in the progenitor, due to the overre- 
action of nuclear processes against perturbations. 

Lai & Goldreich (2000) found an instability in the outer 
supersonic envelope during the collapsing phase. They per- 
formed both analytical and numerical irrotational perturba- 
tion analysis using the conventional polytropic EoS for col- 
lapsing solutions including EWCS of Shu (1977) and post- 
collapse solution of Yahil (1983), and found instabilities in 
supersonic regions. In their derivation, the perturbed fiow 
was assumed irrotational (i.e. without vorticity) and thus 
g— modes should have been excluded. In an example, they 
initiated their calculation with a g— mode perturbation; this 
initial condition appears to contradict the constraints of 
their analytical derivations and numerical analysis. After the 
emergence of a rebound shock, several instabilities have been 
suggested. Intense convective motions may be sustained out- 
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side the neutrino sphere (e.g. Herant, Benz & Colgate 1992; 
Herant et al. 1994). Standing accretion shock instabihty 
(hereafter SASI; e.g. Foghzzo 2001) appears around ~ 200 
ms after the core bounce in some simulations (e.g. Blondin 
& Mezzacappa 2006; Blondin et al. 2003). Burrows et al. 
(2006, 2007a, b) proposed that I = 1 g-modes at ~ 500 ms 
after the stellar core bounce may serve as an agent to ex- 
tract the gravitational energy for the kinetic energy of SNe. 
The roles of these instabilities are still hotly debated and a 
successful SN explosion requires further explorations. 

No significant core instabilities were reported for the 
core-collapse phase. Linear stability analyses were performed 
by GW, Lai (2000) and Lai & Goldreich (2000) for certain 
dynamic flows. GW perturbed their homologously collaps- 
ing solutions and concluded that this collapse is stable for 
acoustic p— modes with the g— modes being neutral convec- 
tive modes. Lai (2000) extended this acoustic stability analy- 
sis to solutions of Yahil (1983) and found no unstable modes 
for 7 ~ 4/3 cases in numerical explorations. Lai & Goldreich 
(2000) studied the stability of collapsing core and claimed 
that the core remains stable in subsonic regions while the 
envelope becomes unstable in supersonic regions. All these 
stable core statements relies on the assumption of a con- 
ventional polytropic EoS. Conventional polytropic gas flows 
and perturbations correspond to a constant speciflc entropy 
and make g— modes just neutral convective modes. 

Our main theme is to examine stability properties of 
a stellar core collapse with a variable specific entropy dis- 
tribution in a homologously collapsing model (Lou & Cao 
2008). By numerical explorations, we classify various pertur- 
bation modes including p— modes and g— modes (e.g. Cowl- 
ing 1941), some of which are oscillatory while others grow 
with time in power laws. As the hydrostatic equilibrium is a 
limiting case of our model, we find connection and evolution 
between perturbations in progenitors of hydrostatic case and 
dynamic collapsing stage. The most interesting result is that 
some stable g— modes in the static case become unstable in 
the dynamically collapsing stage. In particular, the stability 
of each mode now becomes sensitive to the self-similar evo- 
lution of specific entropy. This instability occurs during the 
core-collapse phase, neither before the core collapse nor after 
the core bounce. We speculate implications of such core in- 
stabilities during the collapse phase. For example, the I — 1 
unstable g— mode may lead to the kick velocity of a pulsar. 

SN explosions involve a chain of physical processes. The 
spherical symmetry may be destroyed by a series of insta- 
bilities in the progenitor during the entire SN explosion. For 
example, overstable g— modes of Goldreich et al. (1996) and 
Murphy et al. (2004) may provide seed g— mode perturba- 
tions during the core collapse. Perturbations of our model 
may connect to further instabilities after the core bounce 
and the emergence of an outgoing shock. Our model results 
are highly suggestive and can be tested numerically. 

This paper is structured as follows. Section 2 describes 
general polytropic solutions of spherically symmetric homol- 
ogous core collapse, as a generalization of GW results. Ho- 
mogeneous ordinary differential equations (ODEs) for three- 
dimensional (3D) general polytropic perturbations are ob- 
tained in Section 3. Numerical results are analyzed in Sec- 
tion 4. We consider several aspects of SNe in Section 5 and 
conclude in Section 6. Mathematical details are summarized 
in Appendices A— C for the convenience of reference. 



2 SPHERICAL HOMOLOGOUS STELLAR 
CORE COLLAPSES 

Before starting the time-dependent 3D general polytropic 
perturbation analysis, we first briefiy summarize homol- 
ogous core collapse solutions with spherical symmetry of 
Lou & Cao (2008), making some notational adjustments in 
our model development for the convenience of comparison 
with GW results. The nonlinear partial differential equa- 
tions (PDEs) for ideal hydrodynamics are conservations of 
momentum and mass, Poisson equation for the gravitational 
field and a general polytropic EoS, viz. 



^ + (u.V)u=-ivP-V$ , (1) 

f^+V- (pu) = , (2) 
= 4ttGp , (3) 
P = Ac(r, t)p^, (^^+u-V^log«: = , (4) 

where u(r, t), P{r, t), p{r, t) and <l>(r, t) are bulk flow ve- 
locity, gas pressure, mass density and gravitational potential 
of the flow system, respectively and G — 6.67 x 10~** cm^/(g 
s'^) is the universal gravitational constant, and 7 = 4/3 is 
the polytropic index for a relativistically hot gas. 



It is known that these ideal nonlinear hydrodynamic 
PDEs are invariant under the time reversal operation, 

t^-t, u -u, p, P P, $ ^ $ . (5) 

This property enables us to use an outflow solution to also 
describe a collapse process, which is very important to un- 
derstand homologous core collapse in terms of these expres- 
sions for the time reversal invariance. 

To generalize the analysis of GW model of spherical 
symmetry, we introduce the following time-dependent spa- 
tial scale factor a(t), 

a[t) = p.{t)-'''[^Y" , (6) 

where time-dependent Pc(t) and constant coefficient Kc are 
the values of p and k at the core centre of a massive progen- 
itor star (thus the subscript c). The dimensional vector ra- 
dius r is scaled to a dimensionless vector radius x = r/a{t). 
Consistently, flow variables of u(r, t), P{r, t), p{r, t) and 
<l>(r-, t) are assumed to take on the following forms of 



u — d{t)Ti , (7) 

p=p^{t)f{x)=^^y^\-'f{x), (8) 

(TvGy 

*^^(^)''^-^(.), (10) 



where «:(r, t) is prescribed in the form of Kcg{x) with a con- 
stant Kc- The central mass density Pc{t) is proportional to 
a~^(t) and varies with t. Substituting expressions (f7)l — (|10|) 
into nonlinear PDEs HJ — (SI under spherical symmetry 
with X — r/a{t), we reduce these nonlinear PDEs to a set of 
coupled nonlinear ODEs. First, PDEs ((T| and Q are auto- 
matically satisfled, the latter of which means that g{x) can 
be of an arbitrary form. For example, g{x) = 1 brings our 
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general polytropic EoS back to the conventional polytropic 
EoS with a constant coefScient k — Kc studied by GW. Mo- 
mentum equation then leads to 



a a = — 

X 



P dx ' 3dx 



where the left-hand side (LHS) depends only on t while the 
right-hand side (RHS) depends only on x. For consistency, 
we therefore need to set both sides equal to a constant 4 A/3 
and obtain two separate nonlinear ODEs, viz. 



TVG 
k3 



1/2 



7> > 



(12) 



ODE (|12|) indicates that the spatial scale factor a{t) is either 
a constant independent of t with A = or a power law of t 
being proportional to t^^^. Substituting equation (|13p into 
Poisson equation Q, we derive an ODE for f{x), viz. 



J_d_ 

x^ dx 



x^ d . 4^ 

pTx^^^y 



+ 4/'^ = 4A 



(14) 



which gives a profile of mass density p by eq ([8)l at time t 
and radius r because x is the independent similarity variable 
combining t and r together. Once j(x) is known, other vari- 
ables can all be readily derived. The special case of A = 
leads to the limit of general polytropic Lane-Emden equa- 
tiorQ (e.g. Eddington 1926; Chandrasekhar 1939) and A > 
describes outflows or collapses by the time reversal opera- 
tion. For a necessary check, g(x) — 1 in ODE ((14} reduces 
to equation (16) of GW precisely as expected. 

The 'boundary conditions' for second-order nonlinear 
ODE (|14|l are as follows. The radial gradient of pressure 
should vanish at the centre, i.e. 



VPL= 



: 



4/'(0)+g'(0)=0 



(15) 



where the prime ' indicates the first derivative in terms of 
the self-similar independent variable x. In order to obtain a 
physically sensible solution of f{x) related to the mass den- 
sity by algebraic expression ((8)l, we require an outer bound- 
ary Xb which is the smallest value of solutions f{x) — 0. 
The reason is, if no solution is found for /(s) = at a fi- 
nite a; > 0, it means that the system extends to infinity. 
Therefore, the dimensional velocity u = 2r/(3t) will diverge 
towards extremely large radii. Such a divergent fiow velocity 
is unacceptable in realistic astrophysical gas systems. 

Hence, these boundary conditions determine a contin- 
uous range of A values and the maximum acceptable value 
of A, denoted by Xm hereafter, corresponding to a solution 
f{x) where f'{xt) — also at the outer boundary f{xb) — 0. 
By numerical explorations, solution f{x) for A > Xm does 
not go to zero at a finite x but oscillate with decreasing am- 
phtude with increasing x. Different profiles of g{x) will lead 
to different values of Xm. For the special case of g{x) — 1, 
the A range is < A < Aj\/ = 0.00654376 which was first 



^ The standard Lane-Emden equation is derived by presuming a 
constant entropy in a gas sphere under the self-gravity. For this 
entropy to be a function of r, we refer to the resulting equilibrium 
equation as the general polytropic Lane-Emden equation. 



determined by GW and also confirmed by Lou & Cao (2008) 
with a corresponding fc = /(O) = 4.67047. 

The local polytropic sound speed is defined by 



(11) V. ^ [^y\ [gix)p^^r^ 



(16) 



In the theory of stellar oscillations, it is required that the 
sound speed at the centre approaches a finite value (e.g. 
Unno et al. 1979). We impose the same condition for a dy- 
namic collapse. Hence, g'{0) = and therefore /'(O) — 0. 

Once f{x) is obtained, it is straightforward to calculate 
the total enclosed mass of the collapsing core and the ratio 
between the mean mass density and the central mass density 
in terms of f'{xb) and Xb- The total enclosed mass M of the 
collapsing core is given by 



2 4 3 

M = I 4:TTr pdr — —TvXb 
Jo 3 

3 V TvG / pc 



3/2 



^ 9{n)f'{xb) 

Xb 



(17) 



where the ratio between the mean mass density p and the 
central mass density pc is readily identified with 



— = A g{xb)f'{xb) 

pc Xb 



(18) 



The mean mass density p varies with t because rt decreases 
with t for a homologous stellar core collapse. For A = Am, 
we have f'{xb) — and therefore p/pc = Xm- 



3 THREE-DIMENSIONAL PERTURBATIONS 

We now consider 3D general polytropic perturbations to the 
background self-similar hydrodynamic collapse described in 
the previous section in spherical polar coordinates (r, 9, tp). 
Flow variables including small perturbations are assumed to 
bear the following forms of 



5(^)r(^)[i + /3i(a;, e, ^)r(t)] 



(ttG) 



if,)T(t) , 


(19) 




(20) 




(21) 




(22) 



$ = ^ (^) «W"'[V'(^) + M^, e, ^)T{t)] , (22) 

where the subscript 1 indicates associations with perturba- 
tion terms which are small compared to the background dy- 
namic fiow variables and the free-fall timescale tff is itself 
time-dependent and is written in the specific form of 



'■// 



-nGpc 



-1/2 



9Ay/^ 

2 ' * 



for all A 
for A > . 



(23) 
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In expressions (|19|l — (|22p above, r(t) is a time-dependent 
factor for perturbations and is assumed to bear the form of 



T(t) = exp J tj-ydt^ 



= exp 



1/2 



for all A 

for A > , (24) 



where the value of index parameter p indicates either in- 
crease or decrease as well as oscillations of perturbations 
relative to the self-similar dynamic background flow. Substi- 
tuting expressions (0 - ((221) into nonlinear PDEs ©-(111) 
with the standard linearization procedure, we obtain equa- 
tions governing 3D linear perturbations, viz. 



mwi = -— [V - /iV (gf)] - VV'i 



ffi 3wi ■ V/ /V • wi = 



1 /3i - 3.fi I + wi • V5 = 



(25) 

(26) 
(27) 

(28) 



for vector momentum equation, mass conservation, Poisson 
equation, and specific entropy conservation along stream- 
lines, respectively, where the velocity perturbation Vi — pwi 
and m = p[p + (A/2)^''^]. Therefore we have two values 
p = -(A/8)^/^±(A/8 + m)^''^ for a given A and m. Our defi- 
nition of m here has an opposite sign difference as compared 
to GW definition immediately after their equation (27). For 
A/8 + m < 0, we have a complex conjugate pair for p cor- 
responding to oscillations, while for A/8 -I- m > 0, we have 
two real values of p corresponding to different perturbation 
growth rates for upper plus and lower minus signs in p. 

The angular components of t/^i, /i and /3i can be read- 
ily separated out from the above perturbation equations by 
spherical harmonic^ and wi takes the specific form of 



Wi = ( erWr +egWt-^ + f^V-T-^-^ ] Ylm{8,ip) 

do smO dip ' 



(29) 



By this form of velocity perturbation, the radial component 
of vorticity perturbation is zero, while the ee and com- 
ponents of vorticity perturbation do not vanish in general. 
This allows the possible presence of g— mode perturbations 
as well as convective motions and is distinctly different from 
irrotational velocity perturbations of GW and Lai & Gol- 
dreich (2000). In the following analysis, /i, /3i and xpi only 
describe the radial variations of respective perturbation vari- 
ables with the understanding that the relevant angular parts 



^ To avoid notational confusions, we use index m for spherical 
harmonics Yi^{d, ip) to distinguish from eigenvalue parameter 
m = p[p -\- (A/2)^/^]. The spherical harmonics is defined by 



(21 + 1) (I- m)! 



47r {I + m)! 
m 



1/2 



P"'(cos 6*) exp(im(/j) , 
'I, -(/ - 1), . . . , l-l, I , 
where the associate Legendre polynomial P™ (y) is defined as 



(e.g. Gupta 1978). 



involving the spherical harmonics Yim{9, f) have been sepa- 
rated out. Using equation (|28|) and the angular (i.e. trans- 
verse) components of equation (|25p to eliminate /3i and fi 
in the other three equations, one finally arrives at the fol- 
lowing fourth-order system of homogeneous linear ODEs for 
3D general polytropic perturbations, viz. 



_i_d_ 

dx 



,chfn 
dx 



1(1 + 



3f 



mxwt 



3 dg 



16g 



9d^dg 
4 dx dx 

d 



3m dg 

+ —XWt 

Ag dx 



LA 

x^ dx 



X^Wr) 



/ N 3 do , 
-m— (xwt) + —-rVt- = 
dx Ag dx 

1(1 + i)f ^^df 
X dx 



1 / 

-— mxwt 
9 \ 



(30) 



(31) 



(32) 



Because of the background spherical symmetry, 3D general 
polytropic perturbations are degenerate with respect to the 
azimuthal degree m as expected (note that m and m are two 
distinctly different parameters in our notations). 

Regular boundary conditions at both the centre and 
outer boundary are required to keep perturbations physi- 
cally sensible. They are prescribed as follows. 



Tpl OC x' 



■ IWt 



for X 



( ipi OC X , Sgwrdf / dx — mxwt = ipi for x = xt . 

Boundary conditions at a; ^ 0^ are imposed in order to 
avoid singularity in perturbation solutions at the centre. 
The other boundary condition at the moving radius of a 
collapsing core is to require a zero Lagrangian pressure per- 
turbation there, i.e. dPi/dt + (ui • V)Po = 0. 

For the special case of g{x) = 1 corresponding to a con- 
ventional polytropic gas, these perturbation equations au- 
tomatically reduce to those of GW as expected. Note that 
our Euler equation is written in a vector form because the 
curl of the velocity perturbation field will not vanish for 
a general polytropic gas and thus the stream function ap- 
proach of GW (whose gradient represents the velocity per- 
turbation field) is not sufficiently inclusive especially in view 
of possible g— mode oscillations and convective instabilities. 
However, for a conventional polytropic gas of constant k or 
g{x) = 1, a stream function can be defined and thus GW 
expressed Euler equations in a scaler form by using such a 
stream function. They perturb the stream function instead 
of the velocity field directly. More specifically, if the per- 
turbation stream function takes the form of ^'lyim; then 
wi — V^'i gives 



XWt 



*i 



d^i 

dx 



This illustrates GW result being a special subcase of our 
more general polytropic model description. In other words, 
GW consider only perturbed potential fiows without vortic- 
ity perturbations; this approach suffices for purely acoustic 
oscillations. Likewise, the irrotational perturbation flows of 
Lai & Goldreich (2000) should retain acoustic p— modes but 
exclude gravity g— modes and convective motions. In our 
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perturbation approach, vorticity perturbations are present 
and all possible oscillations for p— modes, f— modes, and 
g— modes are included in the model consideration. 

Several solution properties of this eigenvalue perturba- 
tion problem can be demonstrated. For example, the eigen- 
functions of different eigenvalues m are mutually orthogo- 
nal (see Appendix A). In Appendix B, the eigenvalues and 
eigenfunctions can also be written in terms of the varia- 
tional principle (e.g. Chandrasekhar 1964). In particular, we 
demonstrate in Appendix C that the total energy criterion 
of Chandrasekhar (1939) is not sufficient to guarantee the 
stability of an equilibrium configuration in view of the pos- 
sible onset of convective instabilities for a variable specific 
entropy distribution. 



4 RESULTS OF PERTURBATION ANALYSIS 

In this section, 3D perturbation solutions are divided into 
several classes analogous to the classification schemes of 
Cowling (1941), Cox (1976) and Unnoet al. (1979) for global 
stellar oscillations of static spherical stars. 



4.1 A General Consideration 

Non-radial oscillation modes of a static spherical star have 
been separated into different branches according to their 
respective dominating restoring forces and their frequency 
ranges (e.g. Cowling 1941). 

The p— modes correspond to acoustic oscillation modes 
in which pressure force is the major restoring force. Grav- 
ity modifies such trapped sound waves in several ways. 
The characteristics of such p— modes are: (1) The peaks 
of perturbation functions tend to concentrate towards the 
outer envelope with increasing degree I (i.e. larger / values). 
(2) Their mode frequencies are relatively high, compared 
with eigenfrequencies of other modes such as f— modes and 
g— modes. (3) The more the number of radial nodes in eigen- 
functions, the higher the p— mode eigenfrequencies and the 
more longitudinal the oscillations are. The radial compo- 
nents of perturbed velocity dominate the oscillation in high- 
degree (i.e. I S> 1) modes. 

Another type of oscillatory modes is the so-called inter- 
nal g— modes in which gravitational restoring force takes the 
dominant role. In contrast to p— modes, the g— mode charac- 
teristics are: (1) The maxima of perturbation eigenfunctions 
bury deeply in the stellar interior. (2) Their frequencies are 
relatively low. (3) The frequency goes lower with the in- 
crease of the number of radial nodes (the well-known anti- 
Sturmian property, e.g. Lou 1995). As a limiting case, the 
frequency will approach zero and the perturbation becomes 
nearly horizontal. 

Between p— modes and g— modes, there exist the tran- 
sitional f— modes which have no nodes in both the mass den- 
sity perturbation and the radial component of velocity per- 
turbation. They are essentially surface modes in that per- 
turbations have evanescent behaviours beneath the surface 
layer. When the perturbation degree I becomes very large, 
the perturbation concentrates around the surface layer. This 
mode is closely related to the so-called Lamb waves (Lamb 
1932; Lou 1990, 1991) which propagate in the horizontal di- 



rection and vanish in the vertical direction. These f— modes 
may also be regarded as the lowest-order p— modes 

More specifically, g— modes can be further divided into 
two kinds, namely g^— modes and g~— modes (reading 
g— plus modes and g— minus modes, respectively), accord- 
ing to whether the eigenvalue m is less or greater than zero. 
For global 3D perturbations in static stars, the former class 
of modes is stable while the latter class is unstable; the sta- 
bility property of these modes (especially the g"*"— modes) 
is modified in background dynamic collapses as discussed 
presently. The existence of such two classes of g— modes de- 
pend on the square of the so-called Brunt- Vaisala buoyancy 
frequency M, defined explicitly by 



dlnp IdlnP 



dr 



7 dr 



(33) 



where Q is the magnitude of the local gravitational acceler- 
ation and 7 is the polytropic index of perturbations. 

If is positive everywhere, eigenvalues of m are al- 
ways negative and only g^— modes occur. If is negative 
everywhere, eigenvalues of m are always positive and only 
unstable g~ —modes occur. If M"^ is positive in a certain part 
of the system and negative in another part, both types of 
g— modes may occur (e.g. Lebovitz 1965a, b, 1966 and Cox 
1976 for stability properties of oscillations in static stars). 
Note that the gravitational acceleration always points to- 
wards the centre of the gas sphere (i.e. 5 > by our nota- 
tion). Therefore, the sign oiAf^ is determined by the expres- 
sion in the parenthesis of definition (|33|l . In fact, if the sign 
of the expression in the parenthesis is negative, this region 
satisfies the Schwarzschild criterion for convective instabil- 
ity (e.g. Lebovitz 1965a). Hence, the g~— modes represent 
convectively unstable modes. In the simplest interpretation, 
J\f is the buoyancy frequency associated with a perturbed 
parcel of fluid in a convectively stable medium. Moreover, 
Scuflaire (1974) concluded for oscillations in static stars that 
eigenfunctions of g~— modes can be oscillatory only in con- 
vectively unstable region (i.e. J^^ < 0). On the other hand, 
eigenfunctions of g^— modes can be oscillatory only in ra- 
diative region (i.e. J\f^ > 0). We found similar features for 
perturbations in a homologously collapsing stellar core by 
our extensive numerical explorations. 

Cowling (1941) showed and we readily confirm by defi- 
nition (|33p that for a conventional polytropic EoS P — up'' 
with K being constant for both dynamic core collapse and 
perturbations, the g— modes are simply neutrally stable con- 
vective modes with M'^ = as noted by GW. For this reason, 
GW investigated stability properties of acoustic p— modes 
and suppressed vorticity perturbations in their model anal- 
ysis by introducing a stream function for velocity pertur- 
bation. We note in passing that the situation with two dif- 
ferent polytropic indices for core collapse and perturbations 
deserves a further investigation. 

In contrast, our general polytropic EoS allows fairly 
free options of specific entropy evolution along streamlines 
so that different g{x) gives different profiles of . Conse- 
quently, g— modes can occur in a homologously collapsing 
stellar core and can modify the convective instability crite- 
rion of such a dynamic core collapse in a non-trivial manner. 
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4.2 Numerical Explorations 

The numerical schemes we use to solve this perturbation 
problem are described below. First, given a proper value of A 
parameter and a prescribed g{x), we use an explicit fourth- 
order Runge-Kutta scheme to numerically solve nonlinear 
ODE (|14p for f{x) and then determine all relevant back- 
ground self-similar dynamic flow variable profiles, i.e. tp{x) 
by first-order ODE JTSj, and then u(r, t), P{r, t), p{r, t) 
and $(r, t) by equations (ff)) — (|fOp correspondingly. Hav- 
ing done this, we discretize self-similar perturbation ODEs 
H30p — (|32p by using a proper mesh in order to cast this eigen- 
value problem of ODEs into a matrix eigenvalue problem. 
Because small eigenvalues are practically important, inverse 
iteration method is employed to compute eigenvalues of m 
accurately and to determine the corresponding eigenfunc- 
tions (e.g. Wilkinson 1965). For the purpose of checking, 
we substitute the obtained eigenvalues and eigenfunctions 
into linear ODEs pO|l — (|32|) to compute the residues which 
are sufficiently small. This verification confirms the valid- 
ity of our numerical method. For double check, the Runge- 
Kutta shooting method is also applied (i.e. starting numeri- 
cal integrations near both ends towards the centre) to verify 
the eigenvalues and eigensolutions for perturbations. More 
specifically, this perturbation problem involves undesirable 
diverging solutions towards both ends. We use the numerical 
results obtained from the matrix inverse iteration method at 
two respective points sufficiently close both ends and inte- 
grate towards each other to meet at a mid-point. This avoids 
the numerical difficulty of diverging solutions towards both 
ends at a; = and x = Xb- 

Given a specified profile of g{x) for the specific en- 
tropy distribution with g'{x) = at a; = and a pertur- 
bation degree I value, numerical computations for 3D gen- 
eral polytropic perturbations lead to families of eigenvalue 
m(A) curves. One important feature for these curves as re- 
vealed by our numerical explorations is that none of these 
families of eigenvalue curves intersects with the m = line, 
i.e. each eigenvalue branch m(A) remains either above or 
below m = Q line. Physically, this indicates that when A 
increases from zero to the maximum Am , no oscillations ini- 
tially belonging to g~— mode regime jump into g"*"— mode 
regime across the demarcation m = Q. This is a general 
empirical conclusion on the basis of our very extensive nu- 
merical explorations. 

The physically relevant form of g{x), i.e. the radial 
evolution of specific entropy along streamlines, requires a 
comprehensive understanding of nuclear processes inside the 
high-density stellar core under consideration. In order to ef- 
fectively illustrate essential features of our model through 
numerical explorations, several possible trial distributions of 
specific entropy which have relatively simple analytic forms 
of g{x) are prescribed. It should be emphasized that our 
model analysis regarding g— modes does carry more general 
and important implications for the stellar core collapse. 

In the following, the conventional polytropic results of 
GW model are first carefully examined by our approach and 
then three other types of g{x) are prescribed to explore sta- 
bility properties of p— modes, f— modes and g— modes. 
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Figure 1. The eigenvalues m for the three lowest-order I = 1 
p— modes versus A variation in a conventional polytropic stellar 
core (i.e. constant k or g(x) = 1) in a homologous collapse. The 
upper most one is an artificial mode branch corresponding to a 
displacement of the coordinate origin. The lower two branches 
are acoustic p— modes with zero and one node in the radial veloc- 
ity perturbation, respectively. The eigenvalues of m converge to 
— 25Am/8 in the limit of A — » \m- The dashed curves represent 
GW results for which we have already adjusted the sign difference 
between eigenvalues of ours and those of GW. Our eigenvalues m 
differ from those of GW significantly, even though the correspond- 
ing eigenfunctions agree very well with those of GW. 



0) 

to 

> 
0) 




X 10"' 



Figure 2. The three branches (solid curves) of eigenvalues m 
for the three lowest-order 1 = 2 p— modes versus A variation in a 
conventional polytropic stellar core (i.e. constant k or g{x) = 1) 
in a homologous collapse. They are all eigenvalues for acoustic 
p— modes with zero, one and two nodes in the radial velocity per- 
turbation, respectively and converge to — 25Aj(//8 in the limit of 
A — » Ajvf . The dashed curves represent GW eigenvalues for which 
we have already adjusted the sign difference between eigenvalues 
of ours and GW's. Again, our eigenvalues m differ from those of 
GW significantly, while the corresponding eigenfunctions agree 
very well with those of GW. 
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Figure 3. Acoustic p— mode eigenfunction profiles of our 
x'^f^/'^wt (solid curves) wliicli is equal to r/^/^io (dashed curves) 
of GW for 1 = 2 and A = (i.e. the static conventional Lane- 
Emden sphere) are illustrated here for comparison. The eigen- 
functions of GW are displayed by dashed curves; they are not 
readily distinguishable, because the agreements of theirs and ours 
are extremely close. We suspect that the eigenvalues of GW as 
shown in Figs . [Tl and [2l are in systematic errors. 
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our results in several ways and suspect systematic compu- 
tational errors in the determination of eigenvalue^ by GW. 
Figs. [T]and[2]sliow the first three branches for eigenvalues of 
m versus A for I = 1 and I = 2 p— modes (the branch with 
m > for i = 1 is an exception). For visual comparisons, 
GW eigenvalues are also displayed by dashed curves in Figs 
[1] and (2] We have specifically compared the eigenfunctions 
of ours with those of GW in Figs. [3] and |4] Within numeri- 
cal errors, the corresponding eigenfunctions agree with each 
other very well. We have carefully examined our computa- 
tional procedures and programs and make sure that the first 
three eigenvalues correspond to the first three lowest-order 
p— mode oscillations with 0, 1 and 2 nodes along the radial 
direction. We suspect systematic errors in the computational 
results of GW in terms of their eigenvalues. 

Figs. 13] and [3] show several eigenfunctions for different 
eigenvalues of m with specified A values for the background 
dynamic collapse. In Figs. |3] and [J] dashed curves are those 
of GW results. The coincidence of their curves and ours are 
very good. Fig. [3] shows that similar to classical non-radial 
oscillations of a static spherical star, the nth lowest eigen- 
value of a p— mode correspond to an eigenfunction with n — 1 
number of nodes. Fig.|4]shows shape variations of eigenfunc- 
tions as A varies in succession. The positive m(A) curve in the 
/ = 1 case is an artificial one whose eigenfunction represents 
a movement of the coordinate origin. Another exceptional 
curve noted by GW is for the case of / = 0, representing a 
different choice of the time origin. Except for these two spe- 
cial cases, eigenvalue m(A) curves of all families approach a 
common limit of — 25Am/8 in the limit of A ^ Am. Curves 
of eigenvalue m(A) which are always smaller than — 25Am/8 
must have such a limit. Otherwise, an essential singularity of 
solutions appears for the limiting case of A = Aj\/ (see GW 
for details). Note that eigenfunctions of acoustic p— modes 
become concentrated towards the surface as A approaches 
\m- This singular behaviour of A a/ should be regarded as 
an artifact of the special mathematical character of this lim- 
iting solution. In reality, the Lagrangian pressure perturba- 
tion vanishes at the free surface. A finite sound speed at the 
surface ensures that the normal modes are regular and the 
eigenvalues are isolated 



Figure 4. Here we show I = 2 p— mode eigenfunctions with two 
radial nodes in the transverse velocity perturbation for different 
values of A. These nodes concentrate towards the surface of the 
dynamic stellar core as A approaches \m- The dashed curves show 
GW results taken from their figure 3. Our eigenfunctions agree 
well with those of GW but with a considerable differences in the 
corresponding eigenvalues m as displayed in Fig. |2] 



^.2.1 Comparisons with GW Model Results 

In order to compare with GW results for acoustic pertur- 
bations in a conventional polytropic gas with 7 = 4/3, we 
simply set g{x) = 1 and let all terms containing derivatives 
of g{x) vanish; this simplifies perturbation ODEs (|30|l — 1)32^ 
considerably. Comparing our results with those of GW, the 
eigenvalue curves m(A) for the same oscillation p— modes do 
not coincide with each other well. We have carefully checked 



^ Our model calculations and checks are carried out as fol- 
lows. Using the matrix inverse iteration procedure (e.g. Wilkinson 
1965) which GW also used, we obtain an eigenvalue with its cor- 
responding eigenfunctions. We then insert the eigenvalue and the 
eigenfunctions into nonlinear ODEs I I30I I — I I30I I to verify the re- 
sults. Meanwhile, given an initial value from the calculated eigen- 
functions, we use an explicit fourth-order Runge-Kutta scheme 
to solve ODEs II30I I — II30I I with the calculated eigenvalue again to 
double-check the correctness of the results. In this verification, 
regions near a; = and the outer boundary x = xi, are excluded 
to avoid diverging solutions there. We then calculate and compare 
the counterpart results shown in the figures of GW. After that, we 
find excellent coincidence in eigenfunctions but systematic errors 
in eigenvalues between ours and GW results. 

® As the limiting case of Am only belongs to a purely mathemat- 
ical singularity and makes little sense in physics, in the following 
discussion, we present common features for general values of A 
without any special attention towards this limiting case. 
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Figure 5. Eigenvalues m for various oscillation modes of several 
lowest orders versus the angular spherical harmonic degree I are 
shown here. Three classes of modes, namely p— modes, f— modes 
and g'*'— modes, are identified according to the values of m and 
the structures of their corresponding eigenfunctions. Note that 
the specific entropy g{x) is an increasing function of x in the 
form of expression I I34I I with ei = 0.01 and parameter A = 0.001 
is adopted. The g~— modes do not exist in this case. 



4-2.2 Further Numerical Explorations 

With the physical requirement of g'{0) = 0, we consider 
three types of variable g{x) for numerical explorations. The 
first type is an increasing function of x, the second type is 
a decreasing function of x and the last type has both in- 
creasing and decreasing regions with increasing x. Referring 
to definition (|33[l of Af^ for the buoyancy frequency squared 
which serves as the criterion for the existence of various 
g— modes, we find that the sign of 7V^ closely follows the 
sign of the first derivative of g{x). Thus in order to illus- 
trate g— modes under various situations, such three types 
of g{x) are chosen as representatives. More specifically, the 
three forms of g(x) axe prescribed as even functions of x and 
are given by 



eix 



eix-' + 1 



e2X 



2(62x2 + 1) ■ 

g{x) = 1 + eax^ exp(— a;^/2) 



(34) 

(35) 
(36) 



where £i, £2, and £3 are three parameters calibrating the val- 
ues of respective derivatives and thus the degree of variation 
for the dimensionless g{x). 

Several features of non-radial oscillations are contained 
in the spectrum of eigenvalues m versus spherical harmonic 
degree I for a given value of A and a prescribed form of g{x). 
Figs. [5] E] and [7] schematically present the mode spectra for 
each type of specific entropy g{x) in equations ([34]), (f35)) 
and l|36[l. respectively. Parameters £1 = 0.01, £2 = 0.01 and 
£3 = 0.1 are adopted for our numerical computations. Fig. 
[8] zooms in the neighborhood of the ordinate in Fig.[7l Each 
asterisk * in the four figures represents computed eigenval- 
ues m (i.e. the abscissa) for two or three lowest order modes 
versus the given spherical harmonic degree /. As the back- 
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Figure 6. Eigenvalues m for various oscillation modes of several 
lowest orders versus the angular spherical harmonic degree I are 
shown. Three classes of modes, namely p— modes, f— mode and 
g~— modes are identified according to the values of m and the 
structures of their corresponding eigenfunctions. Only eigenvalues 
m of g~— modes are greater than zero. In this case, the specific 
entropy g(x) is a decreasing function of x in the form of expression 
II35I I with £2 = 0.01 and parameter A = 0.001 is adopted. No 
g+— modes appear in this case. 
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Figure 7. Eigenvalues of m parameter for various oscillation 
modes of several lowest orders versus the angular spherical har- 
monic degree I are shown. Four classes of perturbation modes, i.e. 
p— modes, f— modes, g"*"— modes and g~— modes, are identified ac- 
cording to the values of m. Only eigenvalues of m for convectively 
unstable g~— modes are greater than zero. Here, g{x) takes the 
form of expression I I36I I with £3 = 0.1, neither purely increasing 
nor purely decreasing, and parameter A = 0.001 is adopted. 



ground dynamic core collapse is spherically symmetric, 3D 
perturbations are degenerate with respect to the azimuthal 
degree m (see footnote 2) as expected. Fig. [9] further illus- 
trates the variation of m versus A and reveals that no curve 
of m(A) goes across the line of m = 0. 

By these figures, it is clear to see the variation of 
eigenvalues with changing I. The absolute values of m for 
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Figure 8. The enlargement of region around tlie ordinate axis 
in Fig. [7] It is clear that g+— modes and g^— modes are strictly 
divided by the dashed line m = 0. The g^— and g~— modes 
are distinguished according to the sign of m; the former corre- 
sponds to negative eigenvalues m while the latter corresponds to 
positive eigenvalues m. However, the criterion for convective sta- 
bility is m<— A/8 in a dynamically collapsing stellar core rather 
than m < for a static background. It is therefore possible for 
g+— modes to become convectively unstable. 



Figure 10. As an example, eigenfunctions of a convectively un- 
stable g~— mode for a g{x) in the form of expression II36I I are 
shown. Relevant parameters are A = 0.001 and 1 = 2 and the 
form of g{x) is shown in the figure. The horizontal component 
Wt (dashed curve) is larger than the radial component Wr (solid 
curve). Here, g'{x) is shown by the dash-dotted curve. The peaks 
of eigenfunctions concentrate in the region of g'{x) < corre- 
sponding to Af '^ < 0. 
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Figure 9. Examples of eigenvalues m of ( = 2 p— modes, 
f— modes, g+— mode and g~— mode are shown in this case of 
g{x) = 1 + O.Olx'^ exp(-x2/2). 



p— modes are largest. The higher the order of a p— mode, 
the larger the absolute value of the corresponding m. More- 
over, we find that m is always smaller than —A/8, which lead 
to an imaginary part in the power index of time-dependent 
term r(i) for perturbations. 

Eigenvalues of g— modes, in contrast, lie in the regime 
of smaller absolute value m, being further divided into two 
sub-modes according to the signs of Af^. Figs, [i |6]and[7] 
demonstrate the criteria for the existence of g~ —modes and 
g"*"- modes by /V^, the square of the Brunt- Vaisala buoy- 
ancy frequency A/". In Fig. [5] where the form of g{x) keeps 
A/"^ positive, only g"*"- modes exist and they are thus calcu- 
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Figure 11. As another example of illustration, eigenfunctions 
of a stable g+— modes for a g{x) in the form of expression I I36I I 
are shown. Relevant parameters are A = 0.001 and I = 2 and 
the form of g{x) is given in the figure. The horizontal component 
Wt (dashed curve) is larger than the radial component Wr (solid 
curve). The specific entropy g'{x) is shown by dash-dotted curve. 
The peaks of eigenfunctions concentrate in the central region, 
distinctive from g~— modes and p— modes. 



lated; g~ —modes do not exist as confirmed by our numerical 
explorations. In Fig. [6] where Af^ is negative, only unstable 
g~— modes exist and they are thus calculated; g"""- modes 
do not exist as confirmed by our numerical explorations. 
In Fig. [7] where the sign of Af^ varies in the region, both 
types of g— modes exist and they are calculated by exam- 
ples. Fig. |8] clearly shows the difference between g^— modes 
and g^— modes in terms of eigenvalues of m. 
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Figure 12. As the third example of illustration, unstable 
g+— modes are shown with the eigenvalue m exceeding —A/8. Rel- 
evant parameters and results are shown in the figure. It is clear 
that m + A/8 > which corresponds to an unstable mode. The 
dash-dotted curve is g'{x). The concentration of amplitude peaks 
within the region of g'(x) > 0, i.e. Af^ > 0, and the negative value 
of eigenvalue m both show that it is a g+— mode. 



Recalling the time-dependent factor T(t) by defini- 
tions (|23p and (|24[) and the expression of parameter p — 
-(A/S)^/^ ± (m + A/8)^/^, we find that the criterion for 
convective stability changes with A value of a dynamic core 
collapse. As discussed presently, the first term — (A/8)^''^ 
is interpreted as the gas compression effect. Therefore for 
m + A/8 > of two real values of ±(m + A/8)"'^'''^, one 
mode decays with time t while the other mode diverges with 
time t. For m -|- A/8 < 0, both modes oscillate with time t. 
Consequently, the criterion for convective stability should 
be m < —A/8 rather than simply m < 0. The most inter- 
esting physical consequence is the appearance of unstable 
g"*"- modes of sufficiently high orders. Figs. 1101 1111 and [12] 
demonstrate concrete examples of unstable (j^— modes, sta- 
ble g"*"- modes and unstable g^— modes, respectively. The 
form of specific entropy g{x) is described by expression (|36|) 
for £3 = 0.1. 

The unique f— modes (i.e. Lamb modes), shown in these 
figures and treated as the lowest p— modes, exist only for I > 
2. They separate p— modes from g— modes in the sense that 
the eigenvalues m fall between the eigenvalues of p— modes 
and those of g— modes. An important feature is that the 
eigenfunctions of f— modes contain no nodes in the radial 
component of velocity perturbation. 

Non-radial oscillations of p— modes, f— modes and 
g— modes exist for Z > 1 (in particular f— modes requires 
I > 2). Oscillations for I = modes are purely radial acous- 
tic oscillations as loo = l/(27r)^^^. Since the horizontal com- 
ponent vanishes, one can infer from the primary direction 
of perturbation of different modes that only p— modes exist. 
Our numerical results confirm this inference for a self-similar 
dynamic core collapse in that only stable p— modes exist for 
purely radial oscillations. 

Moreover, several characteristic features of g— modes 
are displayed in Figs. 1101 1111 and 1121 These perturbations 
are primarily horizontal as the horizontal component of per- 



turbation velocity is usually larger than the radial compo- 
nent of perturbation velocity. The peaks of g^— modes are 
trapped deeply in the stellar interior of the collapsing stellar 
core while peaks of g~— mode eigenfunctions appear in the 
region of Af^ < 0. Both have much difference from acous- 
tic p— modes in which the radial component of perturbation 
dominates and oscillation peaks of eigenfunctions appear to- 
wards the surface layer of a collapsing stellar core. 



5 MODEL CONSIDERATIONS FOR SNE 

Homologous collapses have been invoked here to model the 
dynamic phase of stellar core collapse prior to the emergence 
of a rebound shock and the subsequent SN explosion. In nu- 
merical simulations (e.g. Bruenn 1985) and theoretical stud- 
ies (GW; Lai 2000), this phase was claimed to be stable. We 
perform general polytropic 3D perturbation analysis on such 
a homologous core collapse. Analogous to stellar oscillations 
of a static star, 3D perturbations during this dynamic core 
collapse can be classified into several distinct modes. The 
most interesting and important revelation is that convec- 
tively unstable g— modes appear under generic conditions. 
When this happens, several physical consequences for the 
remnant core follow, which may alter the scenario for the 
breakup of spherical symmetry during the core collapse. 



5.1 Radial Evolution of Specific Entropy 

In our model derivations and computations, we start from 
assumptions of a general polytropic gas and the self-similar 
collapse of such a gas sphere under self-gravity. For a ho- 
mologous core collapse, we find that the general polytropic 
EoS is automatically satisfied. A radial evolution of specific 
entropy g{x) needs to be prescribed for a complete solution 
of a homologous core collapse. Moreover, the instability of 
such a dynamic core collapse depends on the radial evolu- 
tion of specific entropy. A brief review for the radial entropy 
distribution in the core collapse phase deems necessary. 

According to Bethe et al. (1979), entropy per nucleon 
(in unit of Boltzmann's constant ks) is ~ 1 — 1.5 at electron 
capture where the mean electron number per baryon Fe is 
~ 0.31. This entropy value varies little such that a conven- 
tional polytropic EoS was adopted by GW. Nevertheless, 
Bethe et al. (1979) noted that multiple nuclear processes 
indeed change this value from this range and this entropy 
variation serves as the conceptual basis of our model. For ex- 
ample, when the breakup of ^®Fe occurs at a mass density of 
~ 5.9 X 10^^ g cm~^, entropy per nucleon varies from ~ 0.5 
to ~ 1.5 at a temperature of ~ 1 — 2 MeV. While a constant 
specific entropy may be an expedient approximation, elec- 
tron capture, neutrino trapping and other physical processes 
certainly make the specific entropy vary with position and 
time. It is this variation of specific entropy distribution that 
leads to g— mode convective instabilities during a stellar core 
collapse before the emergence of a rebound shock. 

During the deleptonization process of a stellar core col- 
lapse, an entropy evolution according to numerical simula- 
tions (e.g. Bruenn 1985; Janka et al. 1995, 1996; Burrows 
et al. 2006) might suggest that specific entropy increases 
with increasing radius, as grossly described by an increasing 
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type of g{x). Bruenn (1985) displayed an increasing distri- 
bution of entropy versus the enclosed mass. As the enclosed 
mass always increases with radius, the specific entropy then 
ascends with increasing radius. More evidence comes from 
simulations for the evolution after the core bounce. Several 
models illustrate the increasing trend of specific entropy dis- 
tribution with radius just after the core bounce (e.g. Janka 
et al. 1995, 1996; Burrows et al. 2006). If the entropy distri- 
bution varies not very much through the core bounce, then 
the specific entropy distribution might also be an increas- 
ing function in the core collapse stage, corresponding to an 
increasing trend of g{x) versus x. 

In short, the analysis of Bethe et al. (1979) indicates 
that entropy can vary in the core collapse phase while other 
numerical simulations (e.g. Bruenn 1985, 1989a, b; Janka & 
Miiller 1995, 1996; Burrows et al. 2006, 2007b) suggested 
an increasing distribution of entropy during the phase of 
rebound shock. Therefore, the overall distribution of entropy 
might be assumed to increase with radius. Having said this, 
the possibility cannot be ruled out that in parts of regions 
the entropy may decrease with radius locally. 

Bruenn (1985) pointed out that entropy decreases with 
density for mass density < 10^^ g cm~^. A turning point 
p* appears here, above which the entropy increases with 
density. The reason is that at a certain nuclear density, the 
energy level of I/5/2 neutron shell in the nuclei becomes filled 
and the electron capture by nuclei is no longer allowed (e.g. 
Fuller et al. 1982). According to this result and the fact that 
mass density is a decreasing function of radius, when p > p* 
at a certain stage, the specific entropy distribution may not 
be a monotonic increasing function of radius. If this analysis 
refiects the physical reality, there are likely some regions 
where the specific entropy decreases with radius. Besides, 
detailed numerical simulations (e.g. Janka et al. 2007) with 
more sophisticated EoS (e.g. Hillebrandt et al. 1984; Shen 
et al. 1998) reveal that although the overall tendency is an 
outward increase, the entropy distribution may have regions 
where the entropy dips slightly. 

The EoS and the specific entropy distribution remain 
open questions. No definitive evidence requires necessarily 
a constant or a monotonically increasing specific entropy 
distribution with increasing radius during the core collapse 
phase. By our computations, variable entropy regions cer- 
tainly lead to unstable g— modes with possibly rapid growth 
rates. These findings bear important implications to pertur- 
bation growths during the core collapse prior to the emer- 
gence of a rebound shock and SN explosions. 



5.2 Definition of a Collapsing Inner Stellar Core 

For the type of SNe involving stellar core collapses, a re- 
bound shock emerges surrounding the centre because the in- 
ner core is drastically compressed and stiffened obeying the 
EoS at nuclear density. After such a core bounce, the central 
neutron-rich core may become a proto-neutron star within 
a mass range of ~ 1 — 3Mq (e.g. Rhoades & Ruffini 1974). 
One expects that in the pre-coUapse stage, there should ex- 
ist a dense central core with a comparable mass collapsing 
inwards to form such a proto-neutron star. Some numerical 
simulations (e.g. Woosley, Heger & Weaver 2002) give the 
central core masses for progenitors of various masses with 



different metallicities. The core mass appears in the range 
of ~ 1.2 - 1.8Mq {Mq = 2 X 10^^ g is the solar mass). 

For a sudden core collapse, the outer layers of the pro- 
genitor may not move in immediately (Burrows 2000). As 
the stellar core rapidly contracts inwards, it may temporar- 
ily detach from outer shells and evolve independently. In 
this sense, we conceive a collapsing core under the self- 
gravity. Meanwhile, a quantitative definition would be desir- 
able. GW used the invariance of the inner core mass to con- 
sider the maximum central pressure reduction percentage 
from a A = configuration to a homologous core collapse. 
However, GW value of ~ 3% for the central pressure reduc- 
tion is much less than the result of ~ 26% in the simulation 
of Bethe et al. (1979). GW defined a core mass as 

= 1.0449 ( — ) Mo , 

where Mo and Kco are the enclosed mass and the value of 
hie for the A = initial static Lane-Emden core while the 
coeflicient 1.0449 comes from the variation of x\pl pc as A 
varies from to the maximum value \m- GW also noted 
that the inner core mass in the computation of Van Riper 
(1978) is ~ 30% larger than that of their definition. 

In reference to GW, we define the following inner core 
mass inside the iron core of the progenitor. Given a form of 
specific entropy evolution gix), we can determine Am and 
then obtain by what factor ci, the value of x\pl p^ increases 
as A increases from to A a/. As this inner enclosed mass M 
is proportional to x\{p/ pc)i^V^ , if we know the value of Kc 
for the two cases of A = and A = Am, the inner core mass 
is then defined by 

/ \ 3/2 
M,, = ci — Mo . 

\K.cO ) 

This definitioijfl has an advantage that it allows arbitrary 
pressure reduction which triggers the SN explosion, though 
the inner core may be very small if the pressure reduction is 
great, especially in the GW cases of conventional poly tropic 
EoS. In our general polytropic EoS characterized by a g(x\ 
the inner core mass can be larger according to this definition. 
In our model development, we actually use parameter 
which is the proportional coefficient between the pressure 
P and p*/'^ at the centre of a collapsing core to estimate the 
inner core mass. According to the results shown in figure 4 
of Hillebrandt et al. (1984) which gave P versus p relations 
for different entropies and values of the central entropy for 
progenitors with different metallicities (e.g. Woosley et al. 
2002), we infer parameter Kc to be > 10^* cgs units. The to- 
tal enclosed core mass is sensitive to the value of Kc because 
of the power- law dependence of hi'J'^ . A form of specific en- 
tropy distribution gixi) = 1-1- O.OOOlx^ exp 

(-iVS) which 

deviates slightly from a constant is chosen here as a demon- 
stration. For Kc = 3 X 10^^ cgs units, the enclosed core mass 
does not exceed ~ O.7M0; for Kc = 5 x lO^'' cgs units, the 



^ In fact, there is yet another definition for the inner core by Yahil 
(1983) in which the edge of the inner core lies on the radius of 
maximum infall velocity. It is not applicable here because Yahil's 
solution extends to infinity and contains both the inner core and 
outer envelope while our solution is valid up to a moving radius 
of zero mass density. 
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enclosed mass is around ~ 1.5Mq; and for Kc = 7 x 10^ cgs 
units, tlie enclosed core mass can reach ~ 2.5Mq. As ex- 
pected, the central part gives the main contribution to the 
enclosed core mass, i.e. materials are highly concentrated 
around the centre. Typically, the chosen form of g{x) does 
not change the enclosed core mass very much unless it de- 
viates from a constant significantly. 



5.3 Properties of Various Perturbation Modes 

5.3.1 Stable oscillations of p— modes and f— modes 

The high-frequency acoustic oscillations of p— modes and 
f— modes may occur in any prescribed form of g{x) and are 
stable during the homologous core collapse. This conclusion 
here is more general than that of GW for a couvcutional 
polytropic stellar core collapse. In our model analysis and 
for convenience, f— modes may be regarded as the lowest- 
order p— modes. Hereafter, we do not distinguish the lowest 
order p— modes and f— modes unless f— modes are necessar- 
ily emphasized. As eigenvalues m of p— modes docroaso with 
increasing orders, if the lowest order mode is stable, then all 
acoustic p— modes remain stable. The local analysis in the 
Wentzel-Kramers-Brillouin- Jeffreys (WKB,J) approximation 
suggests the stability of high-order p— modes (GW). 

A notable feature for perturbation modes in a dynamic 
core collapse reveals that the time-dependent factor does not 
take on a Fourier harmonic form exp(iwt) with uj being the 
angular frequency of a perturbation. Instead, the temporal 
factor acquires a power-law form, except for the limiting case 
of A ^ which consistently reduces to the description of 
harmonic oscillations in a static general polytropic sphere. 
For stable acoustic p— modes, all eigenvalues m found are 
negative and make A/8 + m < where A features the collaps- 
ing core. Therefore, the power-law index parameter T(t) has 
both real and imaginary parts. The real part has a specific 
value of —1/6, corresponding to the adiabatic amplification 
of acoustic waves due to compression of gas collapsing to- 
wards the centre as noted by GW. The imaginary part leads 
to a form of exp(iClnt) with ( = ±(-2m/A - 1/4)^''V3 be- 
ing a real number; this represents an oscillation, but not in 
a conventional sinusoidal form of harmonics. 

For these oscillations of p— modes and f— modes during a 
dynamic core collapse, analyses of GW and ours indicate sta- 
bility. If one includes effects of radiative losses and diffusive 
processes, then some of these acoustic oscillations may be- 
come overstable, i.e. oscillatory growths in a dynamic back- 
ground. Such overstable acoustic oscillations may serve as 
seed acoustic fluctuations for the SASI (e.g. Foglizzo 2001) 
to operate during the subsequent emergence of an outgoing 
rebound shock (e.g. Lou & Wang 2006, 2007). 



5.3.2 Perturbations of g-modes and instabilities 

The most important results of our model analysis are that 
some g— modes which were used to be regarded as convec- 
tively neutral modes under the conventional polytropic ap- 
proximation (GW) can become unstable in a dynamic back- 
ground with a general polytropic core collapse under gravity. 
The occurrence of g"*"— and g~— modes depends on the sign 
of as highlighted in Section 4.1. 



In addition to the power-law factor t^^^^ for oscilla- 
tion amplitudes by the core collapse compression, by us- 
ing the modified onset criterion of convective instabilities, 
i.e. m > —A/8 (see Section 4.2.2), we find that g~— modes 
and sufficiently high-order g^— modes are unstable during 
the core collapse. Here, g~ —modes whose eigenvalues m by 
definition are always greater than zero are unstable. Table 
1 provides information of several illustrative examples of 
g~— modes; eigenvalues m of g"*"— modes approach when 
the order goes higher. As a result, its value must exceed 
—A/8 when the order is high enough. Therefore high-order 
g"*"— modes become unstable. Consequently, unless for a con- 
stant g{x), unstable g— modes always exist, because at least 
one of the two kinds of g— modes occurs for a variable dis- 
tribution of specific entropy. 

It should be emphasized that numerical simulations so 
far have not shown drastic growths of unstable convective 
modes found here. Part of the reason is that there exists 
no systematic study on the stability of the core collapse. 
One would have thought that numerical truncation errors 
and errors in the determination of thermodynamic variables 
from a tabulated EOS are considerable. It would be highly 
desirable to further pursue this problem numerically. 

A physical scenario is advanced below for g"*"- modes. 
When g~— modes do not occur, unstable g"*"— modes may 
give rise to convective instabilities during the core collapse. 
Initially, the inner core of a progenitor remains in a A = 
configuration for which g"'"— modes oscillate stably. When 
an insufficient nuclear energy supply triggers a reduction 
of core pressure, the inner core begins to collapse homolo- 
gously. Sufficiently high order g^— modes become unstable, 
starting to grow and break the spherical symmetry of the 
collapsing core. Such convective instability is limited during 
the core collapse because the growth rate does not exceed 
t~^^^ oc a~^^* because g"*"— modes are defined by m < 0. If 
the spatial scale of the core shrinks by a factor of ~ 10^, 
perturbation grows by a factor of ~ 6. This is comparable 
to the compression effect of a collapsing core. 

We speculate that such g— mode instabilities under fa- 
vorable conditions might break up a collapsing core of high 
density and influence the formation of the central com- 
pact object and its companions, such as binary pulsars (e.g. 
Hulse & Taylor 1975) and planets around neutron stars (e.g. 
Bailes, Lyne & Shemar 1991; Wolszczan & Frail 1992). 

The importance of asymmetry has been emphasized re- 
cently for the core bounce in a progenitor star prior to SN 
explosions. So far, one-dimensional model cannot lead to a 
successful SN explosion. Another consensus is that the neu- 
trino heating mechanism alone also fails to produce a SN 
explosion as the energy of > 10'" ergs appears insufficient 
by one or two orders of magnitude (e.g. Kirauta et al. 2006; 
Burrows et al. 2007a). In recent years, physical mechanisms 
involving two kinds of fluctuations have been proposed to ef- 
fectively extract the available gravitational energy to power 
SN explosions. One is the SASI process, which invokes acous- 
tic fluctuations to transfer energy. The other process relying 
on Z = 1 g— modes appears at several hundred microseconds 
after the core bounce as simulated by Burrows et al. (2006, 
2007a, b). These two energy transfer processes destroy the 
spherical symmetry and make SN explosions anisotropic. 

Regarding the origin of such oscillatory modes or fluc- 
tuations, our proposed instabilities which take place during 
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Table 1. The power-law index -[l/36 + 2m/(9A)]l/2 in the time- 
dopondcnt factor r(t) (see eqns H23| l and l|24p ) for the lowest order 
I = 2 unstable g~— modes varies with A values for the dynamic 
core collapse. There are two values of p resulting from one eigen- 
value m; as the plus sign corresponds to stable modes, we only 
list the index of unstable modes with minus signs in calculating 
index p. In calculations of this example, we find that the decrease 
of eigenvalue m does not exceed 0.6% when A increases from 
(static core) to the maximum value Am = 0.00607. The specific 
entropy distribution g{x) = 1 — O.la;'^ exp(— x^) is adopted ac- 
cording to expression I I36I I. 



A -[l/36 + 2m/(9A)]i/2 



0.001 


-1.663 


0.002 


-1.181 


0.003 


-0.969 


0.004 


-0.843 


0.005 


-0.757 


0.006 


-0.694 


\m = 0.00607 


-0.690 



the core collapse phase before the emergence of a rebound 
shock should have already destroyed the spherical symme- 
try. The condition for the occurrence of such instabilities 
appears generic, only requiring a variable radial specific en- 
tropy distribution. Physically, the g~— modes and unstable 
g"*"- modes lead to convections. By definition (|33[) for the 
Brunt- Vaisala frequency A/", inequality Af^ > is equivalent 
to the Schwarzschild criterion for convective stability in a 
star. The global eigenfunction of a g~— mode describes con- 
vective motions in a homologous core collapse. For a given 
mode, we can calculate the time evolution of the perturba- 
tion. Table [1] gives an example of the power indices for the 
lowest order I = 2 mode under the prescribed form (|36} of 
g{x) — 1 — O.ls^ exp(—2;^). According to Table [2 the power- 
law index varies in a wide range, allowing for both fast and 
slow perturbation growths, indicating the trend that insta- 
bilities go into the nonlinear regime. 



5.4 Comparisons with Stellar Oscillations 

Apparent similarities exist in parallel between stellar oscilla- 
tions of a static star and 3D general polytropic perturbations 
in a homologous core collapse. The classification of different 
oscillatory modes are introduced in reference to stellar oscil- 
lations and the degeneracy with respect to azimuthal degree 
m is expected (see footnote 2). Acoustic p— modes exist for 
all values of perturbation degree I and remain stable. The 
acoustic f— modes exist for I > 2 and also remain stable. The 
criterion for the existence of two different types of g— modes 
remains the same, by using the sign of the square of the 
Brunt- Vaisala buoyancy frequency Af^ defined by expression 
l|33p. Similar to stellar oscillations, amplitudes of g— modes 
concentrate around the core while those of p— modes ap- 
proach the outer layer especially for high-degree I acoustic 
oscillations. More specifically, peak amplitudes of g~ —modes 
appear in the radial regions with AT^ < 0. 

Meanwhile, notable differences from stellar oscillations 
also arise in our perturbation analysis. As the most impor- 
tant features, the onset criterion for convective instability is 



modified in reference to the well-known Schwarzschild crite- 
rion. Consequently, not only g~— modes but also sufficiently 
high order g+— modes become unstable. We provide specific 
examples to demonstrate this novel phenomenon that may 
give rise to several possibilities to inner stellar core collapses 
prior to the emergence of rebound shocks and SN explosions. 
The time-dependent factor no longer takes the exponential 
form but is replaced by a power-law form. 



5.5 Comparisons with Earlier Model Results 

Conceptually, our perturbation analysis during the phase of 
stellar core collapse shows intimate connections to perturba- 
tions before the commencement of stellar core collapse and 
after the onset of core bounce. 

Physically, the origin of perturbations can be fairly nat- 
ural for massive progenitor stars. In terms of stellar evolu- 
tion and dynamics, stellar oscillations of p— modes, f— modes 
and g— modes may well occur in massive stars, such as red 
or blue giants prior to an inner core collapse (e.g. Murphy et 
al. 2004). Nuclear burning in massive stars can also provide 
seed perturbations (e.g. Bazan & Arnett 1998; Meakin & 
Arnett 2006, 2007a, b). Such pre-existing stellar oscillations 
and perturbations associated with nuclear burning serve as 
sources of fiuctuations during the core collapse phase before 
the emergence of a rebound shock around the centre. 

GW studied the acoustic stability of core collapse phase 
for a conventional polytropic gas and noted that g— modes 
are all convectively neutral and p— modes are all stable. Lai 
(2000) extended this polytropic acoustic stability analysis to 
dynamic solutions of Yahil (1983) and concluded that un- 
stable acoustic bar modes (i.e. 1 = 2) exist for 7 < 1.09. 
The instability is caused by an insufficient pressure against 
perturbations due to a soft EoS and may indicate star for- 
mation. A cloud with this instability tends to deform into 
an ellipsoid, in which fragmentation might occur. Except 
for this, no other unstable modes relevant to Yahil solu- 
tions were found in Lai's exploration. Lai also claims that 
Shu (1977) isothermal EWCS is unstable. In the context 
of Type II SNe, he stated that no destructive oscillation 
modes exist in the collapsing core before forming a proto- 
neutron star. These conclusions are based on the assumption 
of a conventional polytropic gas. The fact is that several nu- 
clear processes in SN explosions most likely make specific 
entropy distribution variable. Then g— modes are no longer 
convectively neutral. In particular, their stability now be- 
comes sensitive to the evolution of specific entropy. These 
unstable g— modes should play significant roles for compact 
remnants and SNe. 

Lai (2000) analyzed and claimed the stability of Yahil 
(1983) solution when 7 ~ 4/3; the background dynamic 
flows differ from ours. Lai & Goldreich (2000) found acous- 
tic instability growing during dynamic flows. However, their 
proposed instability only grows in the supersonic region 
where the outer envelope resides (Yahil 1983) while the core 
still collapses subsonically. A main assumption for these re- 
sults is the same conventional polytropic EoS for background 
flow and perturbations. Thus g— modes are convectively neu- 
tral in their analysis and only acoustic p— modes exist. 
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5.6 Several Aspects of SN Explosions 

One can assess consequences of g— mode instabilities during 
the phase of the core collapse prior to SNe. Such g— mode in- 
stabilities can grow and evolve nonlinearly and several plau- 
sible scenarios may be speculated. 

In the presence of a series of core nuclear processes, 
including neutronization and electron capture, fluctuations 
in the specific entropy during the core-collapse phase is in- 
evitable; in particular, such unstable g— modes, i.e. convec- 
tive instabilities destroy the spherical symmetry of a col- 
lapsing core before the emergence of a rebound shock. This 
is a likely mechanism of producing asphericity now seem- 
ingly necessary for SN explosions. Two major consequences 
of such g— mode instabilities are as follows, (i) The formation 
and proper motion of a proto-neutron star can be affected. 
As the collapse becomes aspherical, the proto-neutron star 
may also become aspherical. Under favorable conditions, we 
speculate that such convective instabilities might be violent 
enough to tear a proto-neutron star into pieces, (ii) Such 
instabilities are expected to distort the shape of a rebound 
shock front. Compared with current numerical simulations, 
in which a rebound shock front is initially spherical, then 
becomes non-spherical and bends towards a particular di- 
rection where gas is driven out, our proposed non-spherical 
shock front at the beginning might give rise to considerable 
differences for these numerical simulations. 

There are several conjectures if our proposed instabili- 
ties play a major role in SN explosions. The I = 1 g— mode 
instability may give rise to kicks of proto-neutron stars. The 
I — 2 instability may split the central massive core apart and 
implies a possible formation of binary pulsars; by adjusting 
the parameters in our model, we can have low-order I = 2 
g— mode instability in a central core whose enclosed mass 
is ~ 2Mq. High-order and high-degree g— mode instabilities 
always have a tendency to tear a central core into pieces. 
As a consequence, even though core collapse occurs during 
which neutronization is triggered and energetic neutrinos 
burst out, the outgoing rebound shock may leave behind 
a collection of broken clumps. We suggest this possibility 
because in observations no central compact object is found 
in some SN remnants, including SN1987A (Chevalier 1992; 
Manchester 2007; McCray 2009 private communications). 

One important prediction of our proposed unstable 
g— mode convective instabilities before and during the re- 
bound shock emergence is the resulting turbulent mixing of 
heavier elements in the inner layers of the core for SNe. Ac- 
cording to the evolution theory of massive stars (e.g. Nomoto 
& Hashimoto 1988; Woosley & Weaver 1995), different heavy 
elements (e.g. Fe, Si, O, C etc.) lie in ordered core layers from 
inside out. Without convective instabilities occurring before 
the emergence of a rebound shock, these elemental layers 
should be more or less kept during a SN explosion. In this 
situation, boundaries between elemental shells are expected 
to be identifiable. In contrast, for our proposed g— mode con- 
vective instabilities, the resulting convective turbulence with 
sufficiently fast growths mixes elements Fe, Si, O, C in dif- 
ferent core layers. This convective mixing mechanism of the 
core turbulence should bear observational consequences in 
detecting various heavy nuclear species and mapping their 
spatial distributions in SN remnants. 



5. 6. 1 Speculations on kicks of radio pulsars 

Observational evidence, including high neutron star pecu- 
liar velocities (e.g. Lyne & Lorimer 1994; Cordes et al. 1993; 
Burrows 2000; Arzoumanian et al. 2002; Hobbs et al. 2005), 
the detection of geodetic precession in the binary pulsar PSR 
1913+16 (e.g. Wex et al. 2000), the spin-orbit misalignment 
(e.g. Kaspi et al. 1996), implies "kick" processes by which 
proto-neutron stars gain considerable kinetic energies during 
SN explosions. Three major mechanisms have been pursued 
(e.g. Lai, Chernoff & Cordes 2001). Here, our proposed 1 = 1 
g— mode instabilities in a collapsing core belong to hydrody- 
namically driven "kicks" . We shall not dwell upon the other 
two, viz. electromagnetically driven "kicks" and neutrino- 
magnetic field driven "kicks". 

The dipole modes of Z = 1, including p— modes and both 
types of g— modes, are characterized by possible displace- 
ments of the central core mass as already noted in stellar os- 
cillations (e.g. Christensen-Dalsgaard 1976). As such dipole 
modes involve the core mass motion about the equilibrium 
centre, the core mass gains a certain amount of kinetic en- 
ergy. During the violent rebound shock breakout of a SN, 
this core mass movement enables the remnant compact ob- 
ject to further gain kinetic energy and move away from the 
equilibrium centre. The overall centre of mass should remain 
fixed in space. This is a physically plausible 'kick' process. 
Observationally, proper velocities of nascent neutron stars 
are typically ~ 450 ±90 km s~^ (e.g. Lyne & Lorimer 1994) 
and can reach as high as ~ 1600 km s~^ (e.g. Cordes et al. 
1993; Burrows 2000). For all perturbations with / 7^ 1, the 
boundary condition at j: ^ requires a zero velocity there. 
In these cases, the core mass and the equilibrium centre co- 
incide. Fig. ll4l offers an example of unstable 1 = 1 g^— mode. 
Unstable 1 = 1 g— modes of low orders during the core col- 
lapse phase in a massive progenitor star may give rise to the 
initial 'kick' of a nascent proto-neutron star as joint results 
of nonlinear collapse evolution, rebound shock and SN ex- 
plosion. Here, the requirement of low-order g— modes is to 
enclose a sufficient amount of core mass to be kicked out. 
Our model results indicate that with a plausible distribution 
of specific entropy, low-order I = 1 g— modes can indeed be 
unstable during core collapse. Low-order g— modes ensures a 
sufficient amount of mass inside the inner most node. With 
such instabilities, a sizable mass in the central core can be 
kicked away from the equilibrium centre. 

We may apply this "kick" scenario to a plausible situa- 
tion in which several physical variables are assigned typical 
values for a stellar core. From Hillebrandt et al. (1984) and 
Woosley et al. (2002), the central coefficient Kc is estimated 
as ~ 5 X 10^* cgs unit. The time t is characterized by the 
free-fall timescale ~ {Gp)~^^^ in the order of a few to sev- 
eral tens of millisecond when the mass density p falls in the 
range of ~ 10^" — lO'^* g cm~^. In our model, the time t is 
negative (i.e. time reversal) for a homologously core collapse 
solution; we thus take the initial time to be —50 ms. Here, 
g{x) = 1-1- O.OOlx^ represents a slight increase of specific 
entropy with increasing x. The case of A = 0.006 is close to 
the maximum value Xm of a physically allowed A and also 
leads to a central mass density of ~ 10^" g cm~^ consistent 
with simulation results (e.g. Bruenn 1985). 

Fig. [13] illustrates the initial dynamic core collapse of 
such a situation, including the radius r, the collapsing ve- 
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Figure 13. This example presents the initial dynamic back- 
ground variables for a homologous core collapse such as the radius 
r, the collapsing velocity —u and the mass density p versus the 
enclosed mass M in the unit of the solar mass Mq . 
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Figure 14. We show an example of the I = 1 lowest order 
g+— mode. Parameters for this oscillatory mode are displayed 
in the figure. This is a convectively unstable g+— mode with 
m + A/8 > 0. We can see that the centre of mass has a veloc- 
ity amplitude, which contributes to the kick velocity of the core. 



locity — u and the mass density p versus the enclosed mass 
M. Fig. [13] shows a. I = 1 dipole unstable g+— mode in 
which the central core has maximum perturbation ampli- 
tudes. It is estimated that the sound speed of the central 
core is {■^Pc/ pcY^'^ > 10^ cm s~^ where the subscript c in- 
dicates physical variables at the centre 2; = 0. For a velocity 
perturbation only of a few percent of the sound speed, it 
can still be several hundred kilometers per second which is 
grossly consistent with observed kick velocities of radio pul- 
sars (e.g. Lyne & Lorimer 1994). A somewhat larger velocity 
perturbation amplitude (say 10% of the core sound speed) 
may also give rise to a kick velocity of ~ 1600 km s~^ or even 
higher. We propose that this lower-order I = 1 g— mode in- 



stability is a plausible mechanism leading to kicks of central 
remnant compact objects in the aftermath of SN explosions. 

5.6.2 Speculations on forming binary pulsars 

A few binary pulsars have been detected, including the fa- 
mous PSR1913-I-16 (Hulse & Taylor 1975). Ideas have been 
proposed to explain their formation (e.g. Flannery & van 
den Heuvel 1975). Our I = 2 g— mode instabilities suggest 
an alternative yet plausible origin of such binary pulsars. 
Unstable I = 2 g— modes of low orders developed during 
the stellar core collapse of a massive progenitor might break 
the dense core apart before and during the emergence of a 
rebound shock and eventually give rise to binary compact 
objects. The violent SN explosion drives the two compact 
blobs apart. The spins of the two blobs and the binary or- 
bital motion pick up part of the angular momentum of the 
massive progenitor star. The mass of each components in 
PSR1913-I-16 is estimated to be ~ IMq (e.g. Flannery & 
van den Heuvel 1975). For a central Kc of 7 x 10^* cgs 
unit, the enclosed core mass can be > 2Mq. That would be 
a sufficiently massive core for a possible split into two com- 
pact objects due to nonlinear evolution of unstable I — 2 
g— modes and subsequent core rebounce and SN explosion. 

5.6.3 Speculations on breaking up neutron-rich cores 

We also speculate that under possible and favorable con- 
ditions, the growth of unstable high-order high-degree (e.g. 
^ ^ 3) g— modes and their nonlinear evolution might lead to 
an eventual breakup of a central proto-neutron star after the 
rebound shock emergence and subsequent SN explosion. In 
this scenario, the neutronization does occur during the brief 
core-collapse phase and high-energy neutrinos escape after 
a short moment of trapping, but without forming a coher- 
ent central neutron star in the end. The expected compact 
remnant is actually shredded into pieces by the nonlinear 
evolution of high-order and high-degree unstable convective 
g— modes. SN1987A might be such an example, for which 
no signals from a central compact object have ever been de- 
tected (e.g. Chevalier 1992; Manchester 2007; McCray 2009 
private communications). Unstable high-order high-degree 
g— modes under other conditions might also lead to the for- 
mation of planets orbiting around a neutron star. These 
planets can hardly be those before the SN explosion be- 
cause they can seldom survive under the expansion of the 
progenitor giant and the strong stellar wind as well as the SN 
explosion. These high-order high-degree unstable g— modes 
during the core collapse might help certain amount of gas 
concentrated in isolated blobs, leading to the possible for- 
mation of planets around a new-born neutron star. A few 
neutron stars with planets have been detected observation- 
ally (e.g. Bailes et al. 1991; Wolszczan & Frail 1992). 

We have revealed and emphasized unstable growths of 
g— mode convective instabilities during the inner core col- 
lapse under self-gravity inside a massive progenitor star and 
speculated several physical consequences of such instabilities 
about central compact remnants of SNe. Together with other 
physical processes at different stages, such g— mode convec- 
tive instabilities prior to the emergence of a rebound shock 
can be a necessary ingredient to achieve successful SN explo- 
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sions. Moreover, different modes of oscillations in our clas- 
sification provide sources of perturbations during the emer- 
gence of a rebound shock, for example, the I = 1 g— modes 
discussed by Burrows et al. (2006, 2007a, b) may have well 
originated from our Z = 1 g^— and g~— modes during the 
stellar core collapse phase before the core bounce. Nonlinear 
evolution of diverse initial fluctuations in various stellar core 
collapse conditions is expected to give rise to a diversity of 
possible outcomes for remnant cores of SNe. 



6 SUMMARY AND CONCLUSIONS 

We have systematically investigated physical properties of 
3D perturbations in a hydrodynamic background of self- 
similar core collapse with a general polytropic EoS for a 
relativistically hot gas of 7 = 4/3, as studied in Lou & Cao 
(2008). For both background and perturbations, the two val- 
ues of 7 are taken to be the same; the case of two different 
values of 7 will be explored separately. Analogous to the 
mode classification of stellar oscillations in a non-rotating 
static star (e.g. Cowling 1941), our 3D general polytropic 
perturbations are divided into four distinct classes of modes, 
viz. p— modes, f— modes, g+— modes (i.e. with m < 0) and 
g~— modes (i.e. with m > 0), according to their eigenvalue 
regimes of parameter m in equation (|25|l . 

Stability properties of these different perturbation 
modes are analyzed. Similar to stellar oscillations, acous- 
tic p— modes and f— modes remain stable for 3D general 
polytropic perturbations in homologous stellar core collapse. 
This more general conclusion also confirms the acoustic 
p— mode stabilities claimed by GW although their p— mode 
eigenvalues appear in systematic errors. The temporal am- 
plification factor t~^^^ in the perturbation magnitude is as- 
sociated with the background gas compression during a ho- 
mologous inner core collapse. In contrast, g~— modes and 
sufficiently high-order g"*"— modes are both convectively un- 
stable modes because the onset criterion of convective in- 
stabilities is now shifted from m > for a static general 
polytropic sphere to m > —A/8 for a collapsing general 
polytropic core where A > characterizes the hydrody- 
namic background of a homologously collapsing stellar core. 
The existence of both types of g— modes depends on Af^, 
the square of the Brunt— Vaisala buoyancy frequency (see 
definition I33|l which is determined by the evolution of the 
specific entropy distribution g{x). Meanwhile, above what 
value of perturbation degree I the g— modes are unstable 
also depends on the specific form of g{x). As an example, the 
lowest-order I = 1 unstable g^— mode is shown in Fig. 1141 
The peak amplitudes of g~— modes lie in regions of Af^ < 0, 
which is analogous to those in stellar oscillations. These un- 
stable g— modes lead to growths of convective motions in a 
self-similar collapsing stellar core. Their divergent growths 
scale as power laws in time t (with i < 0) while stable per- 
turbation modes oscillate in the manner of exp(i^lnt). 

In analyzing this perturbation problem, we also realize 
that the global energy criterion of Chandrasekhar (1939) is 
not sufficient to ensure the stability of general polytropic 
equilibria in view of the possible occurrence of convective 
instabilities for variable entropy distributions (Appendix C). 

Compared to possible sources of perturbations proposed 
in earlier models, including the so-called "e— mechanism" be- 



fore the onset of inner core collapse and SASI after the core 
bounce, our g— mode convective instabilities occur during 
the dynamic core collapse. Contrary to earlier theoretical 
notions, the pre-SN stellar core collapse phase is most likely 
convectively unstable due to both types of g— mode insta- 
bilities. This is because an exactly constant specific entropy 
everywhere in a stellar core should be extremely rare in any 
realistic progenitor star (e.g. Bruenn 1985). Therefore, the 
spherical symmetry of a self-similar collapsing stellar core 
should be actually broken up earlier than presumed by most 
previous models of SNe. 

In our scenario, oscillations of the progenitor star serve 
as the most natural source for 3D perturbations. For some 
regions of g'{x) < leading to A/"^ < locally, the unstable 
g~ —modes, of which low-degree modes may have sufficiently 
fast growth rates, will soon dominate and destroy the spher- 
ical symmetry. If g'{x) > everywhere, high-order unstable 
g"*"- modes will grow in a self-similar collapsing core and 
evolve nonlinearly. 

In the presence of inevitable core g— mode convective 
instabilities, several possible consequences may follow. Most 
prominently, the early break-up of spherical symmetry may 
lead to energy concentration in particular directions and 
there is thus no need for a rebound shock to push against 
the entire outer envelope. The low-order unstable I — 1 
g— modes may give rise to the initial kick of a remnant cen- 
tral compact object. Meanwhile unstable g— modes corre- 
spond to the growth of convective motions which may stir 
up heavier elements Fe, Si, O and C in different inner lay- 
ers of the collapsing core inside a massive progenitor via 
convective turbulence so that a mixed distribution of these 
elements to various extents in SN remnants is expected. This 
prediction of our model can be tested by nuclear abundance 
observations in SNe. 

We further suspect that the nonlinear evolution of low- 
order unstable I = 2 g— mode instabilities might disintegrate 
the central proto-neutron star core into two blobs to form bi- 
nary pulsar systems under favorable conditions. Other high- 
order high-degree unstable g— modes might even prevent the 
formation of a coherent central compact object by breaking 
the core into multiple pieces. It is speculated that this might 
happen to supernova SN1987A without signals from a cen- 
tral remnant compact object. High-degree and high-order 
g— mode instabilities may lead to smaller blobs of low masses 
which might eventually form planets around a neutron star 
after a SN explosion. 
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APPENDIX A: ORTHOGONALITY OF 
EIGENSOLUTIONS WITH DIFFERENT M 

We can now prove that eigenfunctions of different eigenval- 
ues m are mutually orthogonal. In the following illustration, 
superscripts and are utilized to distinguish different 
eigenfunctions and eigenvalues. Let us proceed to evaluate 
the integral of 

+ ^ I /(Dy . w'^'dV - 1 fv4'^ ■ w'^'dy . (Al) 
The last term on the RHS of eq. ()A1[) leads to 

' fvi>['^ ■ w^'^dv = I v« V ■ ( / dV 

= ^<^VVf = \ I . v4''dV . (A2) 

The first integration on the RHS of eq (lAlfl can be re- 
expressed as 

lff'(l9A''-^9-^''')^-^'''dV 



Iff fe./P'-V5-wW 



X ( //f ) + 3V/ ■ w^^' ) dV 



(2) 



dv-3 1 /Vf'v/-w(2)dy 



dV 



+ 11 f (Vp.wW)/f)rfF, 



(A3) 



and the second integration on the RHS of eq (|A1|) is equal 
to 

3 / fgfi'^Vf ■ w^'Uv ff['^Vg ■ w^'UV . (A4) 

Note that both / and g depend only on x. Therefore, the 
result indicates that the integral 



is manifestly symmetric in terms of superscripts and 
Consequently, the orthogonality of eigenfunctions is proved. 



APPENDIX B: VARIATIONAL PRINCIPLE 

The variational principle (e.g. Chandrasekhar 1964) can also 
be applied to this eigenvalue problem. Referring to the for- 
mula in the proof of orthogonality in Appendix A, the eigen- 
value can be written as the ratio of two integrals, viz. 



where 

h = I fw^dV 



(B2) 



gf^fUv + l I ivvij'dv 



+ 1 J /'(V/-W) (Vff-w) dV 



+ 2 / r/i(V3-w)dy 



(B3) 



Using the variational principle, we will show that the eigen- 
value parameter m has a stationary property when Ii and 
I2 are evaluated in terms of the true proper solutions. Ac- 
cording to equation (|B1[) . 

(B4) 



5m — {SI2 — mSIi) /Ii , 



where SIi and 5I2 axe the changes in 7i and I2 in response 
to the variation 5w in w. We have 



(B5) 



Sh = 2 I fw ■ SwdV 
and 



5/2 = -2 J gffi5hdV+'^ J vi^^-VS^idV 

+^ f /* (Vg ■ w) SfidV + ^ f f^h {Sg ■ 5w) dV 



+ 1 j fiVf-Sw) (Vg-w) dV 



(B6) 



Keeping in mind of equations (|25p — (|28p . we can write the 
variation 5I2 as 

5l,.2|{-^[V(<;//30 

-.fiV (gf^)] - fV^i} ■ SwdV . (B7) 
By equations (|B5P and (|B7|l . it then follows that Sm = if 

m/V = [V " /iV {gf)] - fVi^i , (B8) 



m — h/h 



(Bl) 



which is precisely equation (|25p for the eigenvalue problem 
of our model formulation. 



APPENDIX C: TOTAL ENERGY OF A 
GENERAL LANE-EMDEN SPHERE 

In the model formulation of Lou & Cao (2008), the situation 
of A = represents the limiting static general polytropic 
sphere. For the conventional polytropic equation of state 
(i.e. a constant specific entropy everywhere independent of 
time t), the governing equation of the static equilibrium re- 
turns to the well-known Lane-Emden equation (e.g. Edding- 
ton 1926; Chandrasekhar 1939). For a general polytropic gas 
sphere with a variable radial distribution of the specific en- 
tropy, we refer to the governing equation with A = as the 
general Lane-Emden equation (see footnote 1 in the main 
text). We emphasize that the linear stability property now 
depends on the radial distribution of the specific entropy. If 
the specific entropy decreases with radius r, corresponding 
to a decreasing g{x), convectively unstable g~— modes will 
develop and drive the static sphere out of equilibrium config- 
uration. The total energy that is equal to the gravitational 
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energy plus the thermal energy together still remains zero 
for 7 = 4/3; however, this situation under the conventional 
polytropic assumption is referred to as being marginally sta- 
ble by Chandrasekhar (1939). 

The hydrostatic equilibrium equations are 

^^4..^, (CI) 

(C2) 

Using equation (|Cip to eliminate Al in equation (|C2|I and 
assuming the radial distribution of P/p~' , one readily arrives 
at the dimensional general Lane-Emden equation whose di- 
mensionless form appears as equation (|14|l with A = 0. 

The total energy of the system can be derived as follows. 
The total gravitational energy is 

Eg = ~ j ^4nr'dr = J ^^Trr'dr (C3) 
and the total thermal energy is 

Et= f , ^ , 47rr^pdr = , , f Pr^dr . (C4) 
J (7-l)P (7-1)7 

The total energy of a general Lane-Emden sphere is 

E = Eg + Et = -{S-y - 4)Et , (C5) 

where finite pressure at the centre and zero pressure at the 
stellar surface are presumed. Consequently, we prove that 
irrespective of the radial distribution of the specific entropy, 
the total energy of a general polytropic Lane-Emden sphere 
remains zero for 7 — 4/3. Evidently, this does not guarantee 
the linear stability of such a equilibrium configuration. 

The simple total energy criterion is not enough to judge 
the stability of the equilibrium configuration. The statement 
of Chandrasekhar (1939) is only valid because the assump- 
tion of a conventional polytropic gas excludes g— modes, 
leaving them as convectively neutral. In fact, linear stabil- 
ity properties rely on the radial distribution of the specific 
entropy. For instance, we find in this paper g~— mode in- 
stability for a static general polytropic Lane-Emden sphere 
with the square of the Brunt- Vaisala frequency Af^ < in 
some parts. We therefore conclude that simply using the 
polytropic index 7 to judge the stability property of a static 
sphere is not enough for a general polytropic Lane-Emden 
sphere. One still needs to perform a detailed linear stability 
analysis for such a hydrostatic equilibrium configuration. 



